Spontaneous creation of non-zero angular momentum modes in tunnel-coupled 

two-dimensional degenerate Bose gases 
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We investigate the dynamics of two tunnel-coupled two-dimensional degenerate Bose gases. The 
reduced dimensionality of the clouds enables us to excite specific angular momentum modes by 
tuning the coupling strength, thereby creating striking patterns in the atom density profile. The 
extreme sensitivity of the system to the coupling and initial phase difference results in a rich variety 
of subsequent dynamics, including vortex production, complex oscillations in relative atom number 
and chiral symmetry breaking due to counter-rotation of the two clouds. 
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I. INTRODUCTION 

A major focus of cold-atom research is coupling mul- 
tiple degenerate Bose gases (DBGs) using atom chips 
and optical lattices fll-[l4|. The results provide step- 
ping stones to future applications, such as interferom- 
eters or processors of quantum information. Since atom 
chips and optical lattices typically generate very high 
(kHz) trapping frequencies, there is growing interest in 
the dynamics of coupled one- and two-dimensional (ID 
and 2D) clouds [15l - [l9| . This interest has also been fu- 
eled by the radically different physics that has been ob- 
served in lower-dimensional systems, such as the suppres- 
sion of equilibration J2(j , quasicondensation [2l[ , and the 
Kosterlitz-Thouless transition [22Tj28j . Since this body 
of work has uncovered such rich dynamics, it is natural 
to wonder how coupled lower-dimensional systems will 
behave, and whether they can reveal a crossover to 3D 
phenomena. Some recent work on ID coupled rings has 
shown that the reduced dimensionality leads to unex- 
pected effects, such as spontaneous population of rotat- 
ing excitations and chiral symmetry breaking [l6j |. Fur- 
ther work is now needed to explore symmetry breaking in 
Josephson junctions, the boundaries of lower-dimensional 
physics, and to establish how double- well interferometers 
will perform in reduced dimensions. 

In this paper, we investigate the dynamics of coupled 
2D disk-shaped DBGs. We find that the reduced di- 
mensionality of the clouds has profound implications for 
the dynamics, because the instabilities of excited states 
observed in three dimensions are suppressed [29|, l3fjj . 
Starting from an irrotational stationary state, and with- 
out introducing any stirring, we observe spontaneous 
occupation of low-lying excitations with non-zero an- 
gular momentum, mediated by the interatomic interac- 
tions [H, [U Hl| • We use a linear stability analysis to 
predict which rotating Bogoliubov modes will be excited, 
and hence identify regimes where we can excite specific 
modes alone by tuning the coupling strength. This tar- 
geted selection of a single rotating mode creates striking 
oscillatory patterns in the atomic density profile. As we 
excite different rotating modes, we uncover a rich pa- 
rameter space. Modes with a high angular momentum 




FIG. 1: (Color online) Schematic constant density surfaces of 
the upper (blue/light grey) and lower (red/dark grey) DBGs. 
Arrows show co-ordinate axes. 



periodically grow and decay exponentially. The growth 
rates of these modes are extremely sensitive to changes in 
the coupling and interaction strengths, and in the relative 
phase of the two DBGs. The growth of modes with lower 
angular momentum becomes unstable due to collisions 
between rotating and non-rotating atoms, disrupting the 
internal structure of the clouds. This leads to a variety 
of subsequent dynamics such as vortex production, os- 
cillations in relative atom number, and chiral symmetry 
breaking due to counter-rotation of the two clouds. 



II. SYSTEM AND METHODOLOGY 

The system consists of two 2D DBGs, referred to as 
the upper and lower DBGs, in the z = +zq and — zq 
planes respectively, as shown in Fig [TJ The two DBGs 
are coupled by a tunnel junction created by the symmet- 
ric double well potential Vd w (z) ■ In the x — y plane, the 
DBGs are contained by the harmonic trapping potential 
V(r) = imw 2 r 2 , where m is the mass of a single atom, 

u) is the trap frequency and r = \J ' x 2 + y 2 . For \z\ ~ Zq, 
Vdw(^) can be approximated by the tight harmonic po- 
tential V£/^(z) — \m\ 2 uj 2 (z =p z ) 2 , where A » 1 to en- 
sure that Xhu > /i. Consequently, atomic motion in the 
z direction is frozen into the single-particle groundstate, 
and the DBG wavefunction becomes 2D [33j. Hence, in 
the weak coupling limit (3^ ]. we may represent the order 
parameter for the two 2D DBGs by the scalar complex 
field 

iP = C(z - z ) X a (p, 0, r) + C(z + z )x?(p, 0, r) (1) 



where the superscripts a and (3 refer to the upper and 
lower DBGs, £ is the normalized single-particle harmonic 

groundstate of (z), 4> is the azimuthal angle, and we 
have introduced the dimensionless time r = tot and the 
dimcnsionless length p = r/a\ 10 , in which a^o = \fh/miu. 
Substituting Eq. ([T]) into the Gross-Pitaevskii equation 
results in two coupled equations for x a ^ (p, 0, r) |35j 
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where the dimensionless quantities 7 = (87rA) ao/flho 
and k = \a\ a f ((z + z ) [d 2 z - 2mV dw (z)/H] ((z - z )dz 
represent the interaction and coupling energy respec- 
tively, and clq is the s-wave scattering length. 

At t = 0, there are an equal number of atoms, No, 
in each well, and X a/0 is the non-rotating groundstate 
of V(r), with chemical potential /iq. For this initial 
configuration, and finite k, there are two possible sta- 
tionary states of Eq. @. These are the ground state, 
defined by x a { r i4>>fy — X^( r i 4>i 0)i with chemical po- 
tential fiQ — \k\, and the excited asymmetric station- 
ary state, henceforth referred to as the 7r-state, defined 
by x Q (r, (/>, 0) = — x^(r,<fi,0), with chemical potential 
A'o + I^l [361 ] . Unsurprisingly, the ground state is stable for 
all coupling strengths. However, the 7r-state shows much 
richer behavior. In three dimensions, the phase discon- 
tinuity may bend, creating vortices via the well known 
snake instability [29|, |3Cj, |371 |38| . This process cannot oc- 
cur in our system, because the reduced dimensionality of 
the disks precludes the movement of the phase discon- 
tinuity. In the following section, we perform a stability 
analysis to identify what excitations may occur. 



III. STABILITY ANALYSIS OF THE tt-STATE 

We perform a stability analysis of the stationary states 
by calculating the excitations of the system using the 
Bogoliubov ansatz |39fl 



x a ^{p,<f>,T)=x a ' p {p,4>^) 
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In the above equation, the m's, v's and mode energies 
Ei have been explicitly labeled with their angular mo- 
mentum quantum number, Z, to derive the four coupled 
Bogoliuvbov-de Gennes equations 
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FIG. 2: (Color) Upper plot: Imaginary component of E Ki i 
(scale shown at top) for angular momentum modes \l\ — to 
4 as a function of Lower surface plots: three examples of 
excited modes with (a) I = 4 and no radial node (\k\ = 1.7), 
(b) I = 1 and no radial node (\k\ = 0.56), (c) I = 1 and a 
single radial node \k\ = 1.35. The vertical height shows the 
atomic density of the mode, and the colors denote its phase 
(blue = 0, green = n, red = 2n). 



where e = 7./V0 and E K> i (= E\ in Eq. [3]) is the energy of 
the Bogoliubov mode for a particular k and I. 

The upper panel in Fig. 2 is a gray-scale plot of the 
imaginary component of E K> i for |2| = to 4 calculated 
as a function of |/c|, for a fixed interaction strength e = 7. 
We are interested in complex eigenvalues because they 
indicate unstable modes which grow exponentially at a 
rate proportional to the imaginary component. The re- 
sults are presented in terms of \l\ because each Bogoli- 
ubov mode has a degenerate pair with equal and oppo- 
site angular momentum. Each shaded region in the figure 
corresponds to a single mode with a particular angular 
momentum and radial excitation. To illustrate this, we 
pick three example couplings of |k = {1.7,0.56,1.35}, 
indicated in the upper panel of Fig. 1 by the vertical 
{dotted, solid, dashed} lines, which pass through shaded 
regions labeled {(a), (b), (c)}. The modes correspond- 
ing to these shaded regions are shown, also labeled { (a) , 
(b), (c)}, in the lower half of Fig. 2 as surface plots, 
whose vertical height denotes the atomic density of the 
mode, and the colors denote its phase. Mode (a) in the 
lower panel of Fig. 2 has four quanta of angular momen- 
tum (I = 4), but no radial node. Mode (b) is similar to 
mode (a) , but has only one quantum of angular momen- 
tum. As we move to larger values of for a given |Z|, we 
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find regions of instability corresponding to modes with a 
greater radial excitation. For example, mode (c) has a 
higher radial excitation than mode (b) because it has an 
additional radial density node (40j . 

The positions of the shaded regions in the upper 
panel of Fig. 2 can be understood by considering the 
limiting case as e — > (no inter-atomic interactions). 
In this limit, the shaded regions shrink to points where 
the coupling strength, k, matches the single-particle 
energies, n + \l\/2, where n denotes the number of 
radial nodes in the corresponding wavefunction. Clearly, 
in this case, excitations can be degenerate. However, 
inter-atomic interactions destroy this degeneracy. In 
general, such interactions shift the shaded regions to 
higher energy. The size of this shift varies, so that there 
are certain values of |k| (for example, \k\ = 1.7, shown 
by the dotted line in the upper panel of Fig. 2) where 
only one mode is unstable, and hence we can tune the 
coupling to excite that mode alone. There are also 
values of k (for example, \k\ — 2.15) where no modes are 
unstable, meaning that the 7r-state is stable. In addition 
to shifting the positions of the shaded regions to higher 
k, increasing the interactions from zero also broadens the 
unstable regions and increases the imaginary component 
of E Ki i within them. This illustrates that, although the 
tunnel coupling provides the energy for the growth of 
unstable modes, the scattering from the 7r-state into the 
excitation is mediated by the inter-atomic interactions. 
We find that e = 7 is a suitable value to provide reason- 
able growth rates of excited modes, without causing the 
unstable regions to become so broad that they overlap, 
so that specific modes can still be accessed individually. 



IV. SIMULATIONS OF THE DYNAMICS 

To simulate the growth of the unstable modes, we 
must model inter-atomic scattering events not cap- 
tured by the Gross-Pitaevskii equation. Hence we add 
non-zero angular momentum fluctuations of the form 

cj a '^r^e~ p2 ^ 2 e il ^ (the single-particle excitations) to 
X a/P [H, then evolve Eq. [2] in time. Here, c\ a/fi) is 
a random complex number such that (|C> |) = and 
(|C (q/ ' 3) |) 2 « lCT 8 0(iVo). These fluctuations are small 
enough that they do not raise the temperature signifi- 
cantly from zero in our model. 

A. Periodic excitation of a rotating mode ([«] = 1.7) 

Figure 3 shows the results of a simulation for |«| = 1.7. 
For this coupling, only a single mode with |/| = 4 [shown 
by surface plot (a) in the lower panel of Fig. 2] is un- 
stable, as indicated by the vertical dotted line in the 
upper panel of Fig. 2. Initially, both clouds have a 
smooth density profile [Fig. 3(a)]. However, as we 



evolve in time, a macroscopic number of atoms are scat- 
tered from the 7r-state via the inter-atomic interactions 
into the / = 4 mode. To preserve angular momen- 
tum, an equal number of atoms transfer into the de- 
generate I — —4 mode, producing the star-like inter- 
ference pattern with eight density antinodes at r « 95 
shown in Fig. 3(b). To quantify the growth of non-zero 
angular momentum modes, we introduce the quantity 
P*/P = l/27r|/x Q//3 (e"^ + e^pdp d(f>\ 2 , which is the 
projection of the upper or lower DBG onto modes with 
angular momentum \l\. We plot log 10 (P") in Fig. 3(d). 
The graph shows that the population of the \l\ = 4 mode 
increases linearly on the logarithmic plot, indicating ex- 
ponential growth, reaching a maximum of approximately 
./Vo/3 at r « 95. At this point, the growth of the I = 4 
mode halts, and atoms begin transferring back into the 
7r-state. The population of the 7r-state does not reach 
zero because the stability analysis presented in Fig. 2 as- 
sumes that all atoms are in the 7r-state, and is therefore 
no longer valid once a macroscopic number of atoms have 
entered an excited mode. Eventually, the DBG returns 
to its original configuration at r « 160 [Fig. 3(c)], at 
which point the process repeats. This periodic transfer 
can be maintained because, due to its large angular mo- 
mentum, the excited mode's peak density is far from the 
origin, meaning that there is little interaction with the 
remaining atoms in the 7r-state. Consequently, the sys- 
tem is able to supresses particle transfer between the two 
wells over a large time, as shown in Fig. 3(e), which plots 
the number imbalance r\ — (N a — N")/2Nq between the 
two wells during the simulation 42] . 

The width of the unstable I = 4 mode in Fig. 2 is 
very narrow, indicating that the formation of the star- 
like interference pattern is a sharply resonant process. 
Hence, the growth rate varies extremely sensitively with 
K, reaching zero if k = 1.7 is changed by only ~ 2%. 
Similarly, the system is also highly sensitive to the rela- 
tive phases of the DBGs. For example, when the above 
simulation is run with an initial relative phase of 0.97T be- 
tween the two clouds, the star-like interference pattern is 
no longer observed. Moreover, an initial relative phase 
of 0.997T causes a 10 % change in the growth rate of the 
star-like pattern. Since this predictable and macroscopic 
pattern is stable, it could be detected experimentally in 
situ using absorption imaging. 



B. Oscillations in relative atom number and vortex 
production (\k\ = 0.56) 

We now tune |«| to 0.56 to excite the lowest possible 
rotating mode, as shown by the vertical solid line in the 
upper panel of Fig. 2. Comparing the surface plots in 
the lower half of Fig. 2, we see that this mode [Fig. 
2(b)] is similar to the one previously excited [Fig. 2(a)], 
except that it has only one quantum of angular momen- 
tum. Consequently, we now expect the formation of an 
interference pattern with two antinodes. Starting from 
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FIG. 3: (Color online) Simulation for \k\ — 1.7. (a)-(c) Atom 
density plot (black = high) of upper (first row) and lower 
(second row) DBG at (a) r = 0, (b) r = 95 and (c) r = 
105. (d) Logarithmic plot of the projection P" of the upper 
DBG onto angular momentum modes with \l\ = 4. (e) Atom 
number imbalance, rj, between upper and lower well. 



the smooth atom density profile at r = [Fig. 4(a)], 
such an interference pattern begins to form, but the two 
antinodes have very different peak densities [Fig. 4(b)]. 
The projection P" shown by the solid curve in Fig. 4(d) 
reveals that, as expected from the stability analysis, the 
population of the |i| = 1 modes initially grows exponen- 
tially. However, as the number of atoms in the \l\ = 1 
modes becomes macroscopic, there is also a significant 
rise in the projection P 2 Q [dashed curve in Fig. 4(d)], 
which is not predicted by the stability analysis. We find 
these unexpected scattering processes when we excite low 
energy modes, which are concentrated close to the cen- 
ter of the trap, and hence collide more frequently with 
the remaining atoms in the 7r-state. Despite these unex- 
pected scattering events, for 10 < r < 40 we may still 
observe one oscillation in the population of the |2| = 1 
mode in Fig. 4(d). However, at r « 50 the oscillation 
breaks down, and vortex pairs appear in the density pro- 
file, as shown in the lower panel of Fig. 4(c), in which 
the position and (anti)clockwise rotation of the vortex is 
indicated by the (minus) plus sign. At the same time, 
we observe the onset of regular, large amplitude (~ 0.8) 
oscillations in the number imbalance rj [plotted in Fig. 
4(e)] due to inter-well transfer. These dynamics reflect 



FIG. 4: (Color online) Simulation for |«| = 0.56. (a)-(c) 
Atom density plot (black = high) of upper (first row) and 
lower (second row) DBG at (a) r = 0, (b) r = 25 and (c) 
r = 75. The positions of vortices are indicated by the black 
(minus) plus, which denotes (anti) clockwise rotation, (d) 
Solid (dashed) line: logarithmic plot of the projection Pf 
(-Bf) °f the upper DBG onto angular momentum modes with 
|i| = 1 (|Z| = 2). (e) Atom number imbalance, rj, between 
upper and lower well. 



that, at this low k, the system has insufficient energy to 
populate other rotating modes and, consequently, it re- 
distributes its energy by performing oscillations in rela- 
tive atom number. This behavior can alternatively be re- 
garded as instability in the symmetric I — groundstate 
mode, which occurs within the shaded region labeled (d) 
in the upper panel of Fig. 2. 

Despite the apparent disruption and the macroscopic 
population of rotating modes in the simulations pre- 
sented hitherto, the system has always conserved zero net 
angular momentum in each DBG, as illustrated, for ex- 
ample, by the presence of a vortex anti- vortex pair in Fig. 
4(c). However the complete system, defined by Eq. ©, 
only requires that the total angular momentum about 
the z axis be conserved, meaning that each DBG may 
break its initial chiral symmetry and obtain a net angu- 
lar momentum, so that = —(if) ^ 0, where (L™^) 
is the angular momentum of the upper/lower DBG. No 
such mode is predicted fom the stability analysis as the 
mechanism for spontanous decay of the 7r-state involves 
two atoms with zero initial angular momentum gaining 
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FIG. 5: (Color online) Simulation for \k\ = 1.35. (a)-(c) 
Atom density plot (black = high) of upper (first row) and 
lower (second row) DBG at (a) r = 0, (b) r = 40 and (c) 
r = 100. (d) Solid (dashed) line: logarithmic plot of the pro- 
jection P" (-Pf) of the upper DBG onto angular momentum 
modes with |/| = 1 (|/| = 2). (e) Atom number imbalance, r\, 
between upper and lower well. Lower plot shows enlargement 
of curve within dashed box in (e) as a (blue) solid curve with 
scale shown left, plus as a (green) dashed/ (black) dotted 
curve with scale shown right. 



equal and opposite angular momentum via an elastic col- 
lision. Since this process is mediated by the non-linear 
term in Eq. (J5]) , it may only occur between atoms in the 
same well. Consequently, we conclude that the break- 
ing of chiral symmetry is a secondary process, caused by 
tunneling of angular momentum between the wells, once 
rotating modes have become populated. We now present 
an example of this behavior. 



Counter-rotation of the two clouds 



1.35) 



We set k to 1.35 (dashed line in Fig. 2) in order to ex- 
cite \l\ = 1 modes with additional radial energy [shown 
in surface plot (c) in the lower panel of Fig. 2] . Figure 



2 reveals that an \l\ = 2 mode is also unstable for this 
coupling, but its imaginary eigenvalue is much smaller 
(lighter gray in Fig. 2) than the corresponding value 
for the |^| = 1 modes. We would therefore expect the 
|/| = 1 modes to grow faster. Our simulation confirms 
this: the initially smooth density profile [Fig. 5(a)] trans- 
forms into an interference pattern [Fig. 5(b)] with two 
antinodes as <f> is rotated through 2-7T, and two antinodes 
in the p direction. The population of this mode [shown 
by the solid curve in Fig. 5(d)] performs one oscillation 
for 15 < t < 45, but by r ks 50 the growing population 
of the |/| = 2 modes [shown by dashed line in Fig. 5(d)] 
has become significant. The occupation of multiple ro- 
tating modes disrupts the DBGs' internal structure, but 
eventually four rough density peaks, encircling the cen- 
tral peak, appear in the atomic density profile, revealing 
macroscopic population of the \l\ = 2 modes [Fig. 5(c)]. 
As this occurs, we observe the onset of oscillations in rela- 
tive atom number, indicated by the periodic fluctuations 
in r] shown in Fig. 5(e). The frequency of the oscillations 
is higher than that shown in Fig. 4(e), because k has been 
increased [34] . Unlike the case presented in section IIVB1 



this is accompanied by a transfer of angular momentum 
between the two DBGs, such that the two DBGs counter- 
rotate. Wc quantify this effect in the enlargement of 
the dashed area in Fig. 5(e), which shows not only r\ 

(solid curve, left axis), but also (L" ) (dashed/dotted 
curve, right axis). The curves illustrate that, although 
— — ), maintaining the net zero angular momen- 
tum of the whole system, (L") and (L@) perform periodic 
oscillations, with a maximum angular momentum differ- 
ence of 1.6/1. Numerically, we find that the lowest energy 
counter-rotating mode with (L^) = —{L") = h has an 
energy of 2No(jiq + 1.3). This explains why macroscopic 
occupation of a counter-rotating mode is not observed in 
the simulation presented in Fig 3, which has an energy 
of 2^0(^0 + 0.56). 



V. CONCLUSION 

In summary, we have shown that the well known decay 
mechanism of a three-dimensional 7r-state [ID, HO, H3, IH| 
is suppressed in tunnel-coupled 2D DBGs. This results 
in a rich variety of physics, such as spontaneous rotation, 
vortex formation, and chiral symmetry breaking due to 
counter-rotation of the two clouds. Our paper builds on 
previous studies of coupled ID rings 16], showing that 
the radial degree of freedom in the 2D disks plays an 
important role. Firstly, excitations in the radial (r) di- 
rection appear in the stability diagram, and may be selec- 
tively populated and experimentally observed. Secondly, 
the clouds may develop a complex internal structure, in- 
cluding vortex anti- vortex pairs. Thirdly, the radial po- 
sition of the peak density, which is fixed in the case of 
the rings, depends on the angular momentum of the 2D 
cloud. As a consequence, disruption arises more readily 
when we populate modes with low angular momentum. 
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Since this work is restricted to zero temperature, it re- 
mains an open question as to how coupling multiple 2D 
DBGs will change their dynamics at finite temperature. 
Finite temperature will introduce a variation of phase 
across the clouds, which might increase disruption, caus- 
ing the counter-rotation to occur more readily. It is also 
possible that coupled systems may provide a link between 
2D Kosterlitz-Thouless behavior and three-dimensional 
Bose-Einstein condensation. 

We note that these systems can be realized with cur- 
rent experimental techniques. Similar set ups have al- 
ready been created by combining a harmonic trap with a 
one-dimensional optical lattice, to investigate quantum 
fluctuation-induced localization Q and the Kosterlitz- 
Thouless transition [22]. Using this method, values of A 
up to 370 have been achieved. A suitable choice of param- 
eters would be, for example, u> = 2ir x 30Hz, A = 10 and 
e = 7. In this case, No ~ 200 87 Rb atoms and the condi- 
tions for two-dimensionality are satisfied. These parame- 
ters are comparable to those in previous experiments Q ■ 



Using optical lattices to create this system also offers 
the intriguing possibility of extendingour work to spinor 
DBGs or arrays of coupled DBGs [JS 0] , which might 
reveal a crossover from 2D behavior to three-dimensional 
band dynamics. 

Finally, we speculate that this system could poten- 
tially be developed as a sensor to detect tiny gravi- 
tational or electromagnetic forces, in a similar manner 
to the production of vortices between merged elongated 
DBGs llHUm. In principle, one of the DBGs could be 
exposed to a tiny force, causing a phase change of 0.017T, 
which could then be detected as a change in the growth 
rate of rotating modes. 
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